Suppression of transverse instabilities of dark solitons and their dispersive shock waves 
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We investigate the impact of nonlocality, owing to diffusive behavior, on transverse instabilities 
of a dark stripe propagating in a defocusing cubic medium. The nonlocal response turns out to 
have a strongly stabilizing effect both in the case of a single soliton input and in the regime where 
dispersive shock waves develop (multi-soliton regime). Such conclusions are supported by the linear 
stability analysis and numerical simulation of the propagation. 



I. INTRODUCTION 

Solitons confined in one dimension (1+1D) are gener- 
ally unstable against periodic modulations in an extra 
dimension (1+2D) [T]. In particular, for the case of defo- 
cusing or repulsive Kerr-like (cubic) nonlinearities where 
the dynamics is governed by the equations of the non- 
linear Schrodinger (NLS) type in two transverse dimen- 
sions, 1+1D solitons are of the dark type [2]. In this case, 
even if 1+1D dark solitons are stable against perturba- 
tions sharing the same dimensionality [3] , the transverse 
instability in 1+2D determines breaking of dark stripes 
into a snake-like shape and subsequent decay into vortex- 
like structures [3H7] , which is the reason why the mecha- 
nism is commonly known as snake instability. Snake-like 
breaking occurs also for bright solitons and other types 
of nonlinear response, and has been the object of recent 
interest [8, 9 . Here we are interested in assessing the 
impact of a nonlocal Kerr response over the transverse 
snake instability of dark stripes. 

A soliton has, by definition, an invariant profile due 
to perfect balance between nonlinearity and diffraction, 
and perturbations grow on top of that over a character- 
istic length scale determined by the inverse growth rate 
of the instability. Potentially, however, snake instabili- 
ties can also affect the propagation of a dark beam in the 
dynamical regime for which the nonlinearity overcom- 
pensates for diffraction. Particularly interesting, in this 
respect, is the strong nonlinear (or weakly dispersive) 
regime, where an initially smooth beam undergoes the 
formation of a dispersive shock wave (DSW), also known 
as undular bore (see, e.g. Refs. [TOHTB] ). In a DSW 
an infinite gradient (catastrophe), which develops owing 
to the initially dominant nonlinearity, is regularized by 
dispersion (diffraction), which determines the onset of 
fast oscillations. These oscillations are basically trains 
of dark solitons [El [16] that fill a region which expands 
beyond the wave-breaking point, i.e. namely a shock 
fan, which is the true signature of the DSW. Such phe- 
nomenon, observed sporadically in the past [T71 [18] , have 



attracted recently a great deal of attention in the area 
of Bose-Einstein condensation (BEC) [T9H2T] and non- 
linear optics [22-25]. In particular experiments in both 
fields demonstrate that DSW can occur from a dark input 
[19] [25]. However, such results suggest that two different 
scenarios can occur, the first of which showing DSW to 
be unaffacted by transverse instabilities [25], in contrast 
with an unstable situation where ultimate decay into vor- 
tices occurs [7l fT9ll2T] . 

In this respect, an important ingredient might be the 
nonlocality of the nonlinear response. In fact, since a 
nonlocal nonlinearity is known to induce suppression of 
modulational instabilities [23[ [29] and transverse insta- 
bilities of bright solitons [30 , the same effect could be 
envisaged for dark inputs. Therefore it is important to 
investigate the onset of snake instabilities of dark soliton 
stripes and their DSW in nonlocal media. In this paper 
we pursue this goal by carrying out the analysis in the 
framework of a sufficiently general differential model con- 
stituted by a NLS coupled to a diffusive equation [23} [26]. 

Concerning the link with previous results, we point 
out that suppression of transverse instabilities of a dark 
soliton stripe has been analyzed in the past, though in 
connection with other mechanisms such as saturation of 
nonlinearities, coupling to other modes or bright solitons 
(vector solitons), or incoherency of the propagating field 
[5] [31]. Viceversa transverse instabilities of DSW are 
nearly unexplored. Investigation of instabilities of shock 
waves have been proposed theoretically in the framework 
of the Burger equation with small viscosity [32]. Con- 
versely, when the dynamics is ruled by weak dispersion 
(versus dissipation), only numerical studies are possible 
since transverse instabilities develop over a strongly dy- 
namic, non-periodic evolution. In this context, our re- 
sults show that an effective suppression of transverse in- 
stabilities of DSW generated by dark stripes is indeed 
possible in 1+2D. 

The paper is organized as follows. In Sec. [II] we briefly 
discuss the model and its solitary waves of the dark type. 
In Sec. |III| we discuss the results of the linear stability 
analysis for a single soliton, while in Sec. |IV| we investi- 
gate the nonlinear stage of the instability showing vortex 
formation over long distances. Finally in Sec. |V| we 
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study numerically the onset of transverse instabilities for 
DSW. 



II. NONLOCAL MODEL AND SOLITARY 
SOLUTIONS 

We investigate how nonlocality affects snake insta- 
bilities in the framework of the following dimensionless 
model, which couples the nonlinear wave equation in the 
paraxial regime for the normalized field u(x, y, z) to a 
diffusive equation which rules the refractive index change 
n n(x,y) 
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a 2 V 2 ± n : 
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(2) 



where a is a free external parameter which measures the 
degree of nonlocality. The local limit is obviously re- 
covered for a = 0, for which Eqs. ([l][2| reduces to the 
standard integrable NLS equation. The model has been 
shown to describe satisfactorily the behavior of thermal 
nonlinearities in liquid solutions [23, 25 , as well as the 
dynamics in plasma [26 j, and has been used as a toy 
model for investigating different aspects of nonlocality 
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phase jump smaller than that of the local soliton. Both 
the modulus and the phase tend in a non-monotonic fash- 
ion to their asymptotic values characteristic of plane wave 
background. In particular for large a they both relax to 
the asymptotic values with damped oscillations. This is 
due to the fact that the eigenvalues of the linearization of 
Eqs. ([!) around the fixed points (plane waves) are con- 
stituted by two complex conjugate pairs, with positive 
(unstable manifold) and negative (stable manifold) real 
part. For fixed <r, solitons exist below a critical maxi- 
mal velocity v max (0 < \v\ < v max ), at which the real 
part of the eigenvalues vanish. We found this upper 
bound to the velocity to be described by the formula 
v max ~ — l)/4o" 2 , which is valid for sufficiently large 
a (a 2 > 1/4). Therefore nonlocal solitons exist for a 
more limited range of velocity compared with the local 
case. Moreover for any given velocity in the domain of 
existence the minimum intensity is larger than the corre- 
sponding local quantity. In other words local solitons are 
darker (their minimum intensity is closer to zero) than 
their nonlocal counterpart with the same velocity. 

Finally, it is important to emphasize that nonlocal dark 
solitons turn out to be stable against perturbations in 
the same coordinate of confinement X, because the nor- 
malized momentum P is always a decreasing function of 
velocity, i.e. dP/dv < [3]. 



A. Dark solitons 

Dark soliton solutions (strictly speaking solitary 
waves) of Eqs. ([l][2| have been investigated in Ref. [28] . 
Here we briefly recall the main results. 1+1D solutions 
are sought for in the form u(x,z) = U(X) exp(z/3Z), 
n(x, Z) = N(X), X = x — vZ being the coordinate in the 
moving frame at velocity v. Here U and N are complex 
and real functions, respectively, which obey the coupled 
equations 



ivd x -p-NjU = 0, 
N-a 2 d 2 x N = \U\ 2 . 



(3) 
(4) 



Without loss of generality we assume the plane-wave 
background of the soliton to have unit amplitude, i.e. 
\U(X = ±oo)| = N(X = ±oo) = 1, which fixes 
the propagation constant to be j3 = —1. For fixed <r, 
with such boundary conditions, solitons are character- 
ized by a single internal parameter, namely their veloc- 
ity v. The soliton family can be constructed by solv- 
ing Eqs. ([3]|4| numerically by means of the relaxation 
method, starting (for any fixed value of velocity v) from 
the explicit solution of the local NLS equation (a = 0) 
U(X) — \f\ — v 2 tanh(Vl — v 2 X)-\-iv, and iterating over 
the parameter a 2 . An example of solution obtained for 
v = 0.6 and a 2 = 5 is reported in Fig. [TJ along with 
the corresponding local solution of the same velocity. As 
shown the soliton has a nonlinear phase profile with a net 
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FIG. 1. (Color online) Example of dark soliton solution ob- 
tained numerically from Eqs. ( 3|4 ) for v = 0.6 and a 2 = 5 
(solid blue line) , compared with the corresponding solution of 
the local case a 1 — (red dotted line): (a) modulus of the 
field; (b) phase of the field; (c) refractive index change N(X) 
across the soliton. 



III. LINEAR STABILITY ANALYSIS 

We proceed to analyze the stability of dark solitons 
against transverse perturbations by using the following 
ansatz 



u = [U(X) + a(X, y, Z)] exp(0Z), 



n = N(X) + b{X,y,Z). 



(5) 
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(a) o 2 = 0, v = (b) a 2 = 0, v = 0.2 (c) a 2 = 0, v = 0.4 
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(d) a 2 = 1 , v = (e) a 2 = 1 , v = 0.2 (f) a 2 = 1 , v = 0.4 
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FIG. 2. (color online) Snake instability gain g computed nu- 
merically vs. transverse wave-number q for different values of 
nonlocality a 2 : (a) stationary case v — 0; (b) v = 0.4. 



FIG. 4. (Color online) Results of numerical simulations based 
on Eqs. (1-2): intensity distribution over the output trans- 
verse plane at propagation distance Z = 10 for different val- 
ues of a and v. Top row, local case a — 0: (a) v = 0; (b) 
v = 0.2; (c) v = 0.4. Bottom row, nonlocal case cr = 1: (d) 
v = 0; (e) v = 0.2; (f) v = 0.4. 
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FIG. 3. (Color online) Overall reduction of the snake insta- 
bility: peak gain g m ax and cut-off (maximal) wave-number q c 
vs. a 2 for three different velocities. 

The perturbations a, b (a <C U , b <C N) obey the lin- 
earized equations 



a - (Na + = 0, 



(6) 

[1 - a 2 (S| + d 2 y )] b = Ua* + U*a. 

Then we write the perturbations in terms of sideband 
pairs with transverse wave-number ±g as 

a = a + {X) ex.p(iqy + iAZ) + a_(X) ex.p(—iqy — i\*Z), 

(7) 

6 = 6 + (X)exp(^2/) + 6_(X)exp(-ig7/), 



where the additional constraint 6_ = b+ must be im- 
posed for the index perturbation to be real. By inserting 
Eqs. ([7]) in Eqs. ([6| and defining the diffusion linear op- 
erator Cd = 1 + <J 2 q 2 — (J 2 d\ and its inverse C^ 1 , we 
obtain from the second of Eqs. ([6| 



6+ = C~ d x (/7a* + U*a) 



(8) 



By substituting this expression into the first of Eqs. (|6|, 
and defining the linear differential operator 



Co 



1 



x 



2 



N(X), 



as well as the unknown array (sideband amplitudes) a = 

[a_|_ a*_] , we obtain the following canonical eigenvalue 
problem 



A& — Aa, 
where the matrix A reads as 



(9) 



A = 



Co - ivd x 



U*C^U* 



uc^u* 



-UC-/U 



-Co - ivd x 



U*Cj 1 U m 



We solve the eigenvalue problem (|9]) numerically by intro- 
ducing a suitable discretization of the operators Co, Cd 
along the X axis. Instability occurs for negative imagi- 
nary eigenvalues, entailing a perturbation that grows ex- 
ponentially with rate g = — Im(A). The eigenfunctions 
(not shown) that correspond to unstable eigenvalues have 
even symmetry, which associated with the odd symme- 
try of the soliton, leads to a snake type of instability. 
The dependence of the snake instability gain g on the 
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wave-number q is displayed in Fig. [2] for different values 
of a both in the stationary (v = 0, still black solitons) 
and non-stationary case (v = 0.4, moving solitons). As 
shown, when nonlocality becomes progressively more ef- 
fective (i.e., when a is increased), a substantial suppres- 
sion of the instability is observed. In fact, both the peak 
gain g max = max[#(g)] and the cut-off wave-number q c 
(below which g ^ 0) decrease as a grows larger. The de- 
pendence of peak gain g m ax and cut-off wave-number q c 
on the nonlocal parameter <j, for three different represen- 
tative velocities v = 0, 0.2, 0.4, is summarized in Fig. [3j 

In order to verify that the snake instability is effec- 
tively suppressed during propagation of a dark stripe, we 
resort to numerical integration of Eqs. (l)-(2) by split- 
step method. We use initial data uq = u(x, y : z = 0) = 
U(x) 2/), £(x,y) denoting henceforth a noise term 

(\£>( x iy)\ ^ l^( x )l) with random (normally distributed) 
modulus and phase. The results of the simulations, com- 
paring the local (a = 0) and the nonlocal (a = 1) case, 
for three different values of velocity are reported in Fig. [4] 
In the local case, one observes the onset of snake instabil- 
ities for any value of velocity (the instability is weaker for 
larger velocities), as clearly shown in Figs. [4Fa-c). How- 
ever, in the nonlocal case, as shown in Fig. Hfd-f), the 
instability is effectively suppressed and the stripe sym- 
metry of the input is preserved. 

It is worth emphasizing that we find a similar stable be- 
havior (results not shown) also for bound states formed 
by two soliton stripes [28]. However, in this case, it is 
improper to talk about the stabilizing role of nonlocal- 
ity, since these states do not survive in the limit of local 
nonlinearities, being peculiar to the nonlocal response. 



IV. VORTEX SOLITONS AND THEIR 
GENERATION 



Having shown that nonlocal dark stripes remain stable 
over characteristic distances where their local counter- 
parts break down, one might wonder about propagation 
over much longer distances (even though these can be 
hardly reached in practical experiments). More specifi- 
cally it is natural to ask whether the weak residual insta- 
bility determines a decay scenario qualitatively similar to 
the local one, where an unstable dark stripe leads to the 
generation of vortex solitons [H|5]. In order to investigate 
this point, we have first searched for vortex solutions of 
the nonlocal model, and then explored whether they are 
generated through the long-range instability of the dark 
stripe. 

Stationary vortex solitons are characterized by a field 
u(r,(p,z) = U(r) ex.p(irrup + i/3z), where the topological 
charge m fixes the azimuthal dependence (helicity of the 
phase) and r = ^ x 1 + y 2 is the radial variable. The 
associated index n(r, <p) = N(r) has only radial depen- 
dence since it depends on intensity. U and N obey the 
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FIG. 5. (Color online) Radial intensity profiles (red solid line) 
of nonlocal vortex solutions and relative index profile (green 
dot-dashed line) for different topological charges m = 1,2,3, 
and cr=l. For comparison we report also the local intensity 
profile (dashed blue line). 



following equations 



1 (d 2 r + -d r -m 2 

2 \ r 



U — NU = 0, (10) 



N-a 2 ^ + -d r jN=\U\\ (11) 

which can be solved numerically by the relaxation 
method. An example of the radial profile of such so- 
lutions, obtained for a fixed value of nonlocality a = 1 
and different helicities m = 1, 2, 3, is displayed in Fig. [5| 
As shown, the nonlocal profiles are qualitatively similar 
to the corresponding local vortex solitons, except for the 
fact the index change does not vanish in r = (where 
the field instead vanishes) because of the averaging due 
to the nonlocality. 

In order to show that these vortex solitons are stable 
modes of the system and appear spontaneously due to the 
residual instability of the dark stripe, we have extended 
the simulations of Fig. [4] to much longer distances. A 
typical result, obtained for a = 1, is reported in Fig. [6] 
As shown in Fig. [6^a), at sufficiently long distances, the 
soliton breaks down into an irregular sequence of dark 
spots. A closer look to two adjacent spots [see insets in 
Figs. [6^b-c)] reveals that they are indeed vortex pairs 
with opposite charge (helicity). The radial profile of the 
spontaneously generated vortex modes, though affected 
by the noise, are closely reminiscent of the solutions of 
Eqs. (10-11), as clearly shown in Fig. [HJd). 

The results shown in Fig. fallows us to conclude that 
the residual instability leads to a scenario dominated by 
spontaneous formation of vortex pairs, similarly to the 
local case. In this case, however, much longer distances 
are needed to observe the phenomenon because of the 
substantial weakening of the instability due to the non- 
locality. 
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FIG. 6. (Color online) Results of numerical simulations based 
on Eqs. (1-2) of long term evolution of an input stationary 
dark soliton stripe, (a) intensity distribution over the output 
transverse plane at Z — 80. (b-c) insets showing a portion 
of the transverse plain including two adjacent vortex solitons: 
(b) intensity distribution; (c) phase distribution, showing op- 
posite helicity. (d) noisy intensity profile along x (blue solid 
line) of the upper vortex in (b) as seen at fixed y = 5, com- 
pared with the theoretical profile (red dashed line) of the vor- 
tex mode as obtained from Eqs. (10-11). Here a — 1. 




FIG. 7. (Color online) DSW ruled by the local 1+2D NLS 
equation, contrasting the ideal noise-less case (a-b-c) with the 
evolution in the presence of noise (d-e-f). (a-d) Shock fan at 
y — 0; (b-e) Output transverse plane; (c-f) Output intensity 
profile (solid line) at y — superimposed to the input (dashed 
line). Here e = 0.05, and the soliton stripe input is U(x) = 
tanh(x). 



INSTABILITIES OF DISPERSIVE SHOCK 
WAVES 



When the intensity becomes much higher than that 
needed to support a soliton, nonlinearity overcompen- 
sates for diffraction, and eventually leads a smooth input 
to develop an infinite gradient (catastrophe) which is reg- 
ularized by diffraction, forming a DSW. In the case of a 
hyperbolic tangent input the catastrophe occurs in the 
beam center, and post-shock oscillations emerge, in the 
dynamics governed by the local 1+1D NLS, in the form 
of a soliton train, with pairs of dark solitons emerging at 
progressively larger angles (transverse velocity) around 
a central immobile soliton. A qualitatively similar sce- 
nario takes place in the nonlocal case [25 j. Here we in- 
vestigate the stability of such DSW structures against 
transverse perturbations in 1+2D. The highly dynami- 
cal non-periodic character of the 1+1D evolution does 
not allow to extend the linear stability analysis carried 
out for single solitons, and hence we resort to numerical 
simulations. In particular, we restrict ourselves to input 
conditions uq = u(x, y, z = 0) = p [U (x) + y)] , where 
U(x) is the 1+1D (either local or nonlocal) soliton, and 
p ^> 1 is the factor which measures ovecompensation (i.e., 
deviation from soliton case). However, in order to make 
evident that, here, we operate in the weakly dispersive 
limit, we introduce the scaling commonly adopted in the 



theory of DSW [25], casting Eqs. (l)-(2) in the semiclas- 
sical form 



is—- + — - nip = 0, 
oz 2 



(12) 
(13) 

which is obtained from Eqs. (l)-(2) with the positions z — 
Z/e, h = ne 2 , ip — eu, where e = 1/p <C 1. Such scaling 
allows us to keep the input fixed to uo = U(x) + £(x,y) 
and increase the nonlinearity by decreasing s. In the fol- 
lowing we present specific examples obtained for e = 0.05, 
using split-step method with a typical transverse grid of 
2048 x 256 points (x x y) over the range of interest. Note 
that, in order to employ periodic boundary conditions 
along the confinement coordinate x, the calculations are 
carried out over a window much larger than that shown in 
the figures below, containing two solitons with opposite 
phase and sufficiently far from each other, so to ensure 
negligible interaction. 

The formation of a DSW occurring in the local limit 
<j = from a black soliton initial shape (U (x) = tanh(x), 
corresponding to v = 0) is illustrated in Fig. [7] For 
reference we display in Fig. |7Fa-c) the ideal DSW ob- 
tained by integrating Eqs. (|12|)-( 13) in the absence of 
noise. We recall that the occurrence of wave-breaking 
can be described in terms of hydrodynamical variables 
h = |^| 2 (water height or density of an equivalent fluid) 
and V — (f x (velocity of the water or the equivalent 
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break [see Fig. (7^e)]. Solitons with larger darkness (i.e., 
dip of intensity closer to zero) break down at shorter dis- 
tances, consistently with the fact that the gain decreases 
for larger velocity, as can be inferred by comparing the 
gain curves in Fig. 2(a) and 2(b). The onset of the insta- 
bility makes the shock fan, once seen on the plane y = 0, 
irregularly filled with erratic filaments which change ve- 
locity in unpredictable manner, as clear from Fig. J7^d). 
The instability not only results into a strongly irregular 
train along x [compare Fig. [7^c) and Fig. [7jf)], but 
also into the fact that such train looks different when 
seen at different elevations (i.e., planes y = constant). 
Importantly, however, in a nonlocal medium, the sta- 
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FIG. 8. (Color online) As in Fig. 5, showing stabilization of 
a DSW against transverse instability, in the nonlocal regime. 
Results from numerical integration of Eqs. ( 12 )-( 13 ) with a = 
1 and e = 0.05. (a) Shock fan at y = 0; (bJOutput (z = 8) 
transverse plane; (c) Output intensity profile (thin solid blue 
line) and refractive index change n(x) (thick solid red line) at 
y = superimposed to the input (dashed line). 
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fluid), which obey the hyperbolic equations (see, e.g. 

Refs. 12a nana GS] ) 



h z + (hu) x =0, 
u z + uu x + h x =0, 



(14) 
(15) 



obtained at leading order b y ap plying the transformation 
ip = Vh exp(i(,2 / e) to Eqs. (I12|-(13) in the local limit. In 
the case shown the gradient catastrophe occurs at the 
point of vanishing intensity x = 0. The catastrophe be- 
comes manifest at finite distance [z = z s ~ 0.75) as a 
vertical front in the velocity V accompanied by a cusp in 
the intensity (density) [25]. Beyond such wave-breaking 
point, diffraction along the x— coordinate regularizes the 
infinite gradient, and a DSW becomes apparent as a fan 
[see Fig. (7^a)], which is progressively filled with trans- 
verse oscillations taking the form of narrow gray soli- 
ton filaments. Asymptotically (z ~ the number 
of emerging solitons is 2p — 1 for p = e~ x integer, and 
the velocities of symmetric pairs tend to accumulate to- 
wards linear velocity \dx/dz\ = 1 [16 . In the absence 
of noise the solitons in the fan emerge onto the trans- 
verse plane as stripes, as shown in Fig. (7^b). However, 
such scenario changes radically when noise is introduced, 
as illustrated by Fig. [7^d)-(e)-(f). Transverse instabili- 
ties rowing from noise cause indeed the soliton stripes to 



FIG. 9. (Color online) DSW dynamics in a local medium 
(a = 0) arising from a gray soliton shape uo = 
y/1 — v 2 tanh(Vl — v 2 x) + iv, with v — 0.5. (a-b) Results 
from numerical integration of Eqs. ( |12| -( 13 ) in 1+1D, with 
e — 0.05: (a) the twin shock fan in the x — z plane. The 
dashed line indicates the natural velocity of the input in the 
limit £ — 1 for which the input itself is a one-soliton; (b) 
Output intensity in the soliton trains (solid line) compared 
with the input intensity (dotted line), (c-d) Comparison of 
the results from Eqs. (12)-(13) (solid line; blue online) with 



those obtained in the hydrodynamic limit from Eqs. ( 14 )-( 15 ) 
(thick dashed line; red online): (c) equivalent fluid density 



h 



(d) equivalent fluid velocity u — cp x vs x. 



bilizing effect of the nonlocality on dark solitons allows 
us to envisage the stable formation of a regular DSW. 
This is demonstrated in Fig. [8j which is obtained by 
integrating Eqs. (12)-(13) with a = 1, e = 0.05, and 



input U(x) corresponding to the exact nonlocal soliton 
shape (quite similar results, not shown for brevity, are ob- 
tained with hyperbolic tangent input). First, as shown in 
Fig. [8^ a), the wave-breaking scenario and the shock fan 
look qualitatively similar to that of the local case without 
noise, except for an increased shock distance z s (already 
pointed out in Ref. [25]), and an adiabatic change of the 
width of solitons in the train. In this case, the details 
of the 1+1D evolutions could be described by means of 
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Whitham modulation theory, whose application does not 
require integr ability of the governing equations. While 
this is clearly beyond the scope of this work, our goal 
here is to emphasize that the output transverse pattern 
of intensity exhibits no transverse break up [see Fig.[8jb)] 
in spite of the noisy input. Conversely, the DSW develop 
in such a way that soliton filaments in the fan remain 
invariant stripes along y. Once seen at a fixed elevation 
y the DSW looks as a soliton train which maintain close 
similarity to the train that one obtains when the propaga- 
tion is strictly 1+1D [25 . In summary, the comparison of 
Fig. [7] and Fig. [8] allows us to conclude that the nonlocal- 
ity can play a crucial role, enabling effective suppression 
of snake instabilities and observation of effective 1+1D 
DSW in a 1+2D setting. The effectiveness of such sup- 
pression mechanism is witnessed also by the fact that we 
observe the same scenario if we increase the initial noise 
level substantially (i.e., one order of magnitude) with re- 
spect to the local case reported in Fig. [7^d-f). Such 
suppression is effective also for a < 1 (results not shown) 
and becomes more and more pronounced as a increases. 
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FIG. 10. (Color online) As in Fig. [9] in 1+2D, showing the 
onset of the instability growing from noise, (a) Shock fan 
at y = 0; (b) Output transverse plane; (c) Output intensity 
profile at y = (solid line) compared with the input (dashed 
line). 

Finally, we have extended the numerical calculations 
in order to investigate the case of a dark input with 
non-zero minimum intensity. First, we point out that, 
in this case, the wave-breaking scenario that leads to the 
appearance of the DSW changes qualitatively, as illus- 



FIG. 11. (Color online) Stabilized DSW in the nonlocal case 
with a — 1. Results from numerical integration of Eqs. (12)- 
(13) with e — 0.05 and input shape characteristic of the non- 
local soliton U(x) corresponding to v = 0.5. 



trated in Fig. [9] in the local case. We employ an in- 
put waveform corresponding to a gray soliton U(x) = 
y/l — v 2 tanh(Vl — v 2 x) + iv of the integrable NLS equa- 
tion. In this case, as shown in Fig. [9^a), a central black 
still soliton no longer appears. Conversely two distinct 
fans develop, which correspond to two non-symmetric 
soliton trains which split up on propagation, having 
net velocities of opposite sign [see Fig. |9jb)]. The two 
fans originates from regularization of two distinct shocks 
which develop trough gradient catastrophys at points 
x = x si located asymmetrically on the right (x s > 0) 
and left (x s < 0) of the input intensity dip (in x = 0), 
respectively. These shocks develop at slightly different 
distances z s , as one can predict from numerical solutions 
of the hydrodynamic limit [Eqs. ([l4|-( 15)] [? ]. Snap- 
shots of the field at distance z s = 1.52 corresponding to 
the formation of the right shock at x s > are shown 
in Fig. [9jc-d), while the left shock occurs at a slightly 
longer distance. In particular Fig. [9^c-d) compare the 
density h = \iJj\ 2 and velocity u = <p x of the equivalent 
fluid, as obtained by integration of reduced hyperbolic 
[Eqs. @-([l5])] and the full model [Eqs. re- 
spectively. As shown, the agreement is remarkably good 
up to the distance considered, beyond which the dynam- 
ics ruled by the full model deviates qualitatively starting 
to exhibit the soliton train features illustrated in Fig.[9^a- 
b) . It is worth pointing out that the asymmetry between 
the right and left shocks is due to the soliton phase (im- 
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plicitly considered in the initial condition U(x)), which 
is responsible for the initial profile of fluid velocity (i.e., 
the derivative of the phase with respect to x) reported in 
Fig. [9jd) as a dotted line. If such phase is suppressed, 
the two shocks turn out to be symmetrically located in 
x around the origin, and develop at the same distance. 

As shown in Fig. [TUJ we find that the 1+1D scenario 
illustrated before changes once the propagation occurs 
in 1+2D, similarly to the case of black soliton input. 
In fact, in 1+2D, the two DSW fans exhibit transverse 
break-up [Fig. [To^b)] due to the onset of the instability 
from noise. As a result, once the dynamics is seen on 
a plane y = constant, the solitons in the fans show er- 



ratic changes of velocity, as can be seen in Fig. 10 ^a). 
Note that the instability develops from the most unsta- 
ble soliton stripes (i.e., the inner ones) in the two fans. 
When seen at given distance z and transverse location y, 
the soliton trains appear irregular in x, as displayed in 
Fig. 10 'c). Moreover the pattern shown in Fig. 10 ^c) 
changes from plane to plane y = constant. 

Importantly, however, when nonlocality becomes effec- 
tive the twin shock fan is stabilized. A typical example, 
corresponding to the internal parameter value v = 0.5, 
already employed in the local case [Fig. [9] and Fig. [To], 
is illustrated in Fig. 11 In this case the input U(x) cor- 
responds to the exact soliton shape obtained with e = 1, 
for the external parameter a = 1 and v = 0.5. In spite 
of a noise level much higher than the corresponding one 
in the local case, the DSW does not break up, showing 
that the suppression of the instability is effective also in 
this case. Similar results are obtained for different val- 
ues of soliton velocity v, allowing us to conclude that the 



suppression of the instability due to nonlocal effect is a 
general phenomenon that permits to observe true 1+1D 
(stripe) DSW dynamics also in 1+2D under different ini- 
tial conditions. 



VI. CONCLUSIONS 

In summary we have shown that nonlocality plays a 
beneficial role in the propagation of dark stripes, allow- 
ing for a substantial suppression of transverse (snake) 
instabilities for both a single soliton and dispersive shock 
waves characterized by the emergence of multiple dark 
soliton stripes. The residual weak instability leads, over 
long distances, to decay into vortex pairs, as numerical 
simulations of the nonlocal model seem to indicate. 

We point out that the nonlocal mechanism of suppres- 
sion investigated in this paper, could cooperate, under 
suitable conditions, with other mechanisms previously 
investigated [31]. Moreover, stabilization due to nonlo- 
cality of dark stripes is expected also with focusing non- 
linearities, when the signs of second derivatives in the 
transverse plane are opposite (hyperbolic case). 
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